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Ordered graphene has been extensively studied. In this paper we undertake a first density functional study of topo- 
logically disordered analogues of graphene, in the form of a random network, consisting predominantly of hexagonal 
rings, but also including pentagons and heptagons. After some preliminaries with crystalline material, we relax various 
random network models and find that the presence of carbon pentagons induce local curvature, thus breaking the initial 
planar symmetry, in some analogy with the case of fullerenes. Using density functional theory to calculate the total 
energy, we find that while the planar state is locally stable, there is a puckered state that has lower energy. The scale of 
the puckering is consistent with that expected with local maxima and minima associated with pentagons surrounded 
by larger rings; forming local "buckyball domes". 
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1 Introduction Graphene is among the hottest top- 
ics in current condensed matter science. A vast amount 
of work on many aspects of crystalline graphene has ap- 
peared. In this paper we take a different tack: we explore 
the role of topological disorder in amorphous graphene. 

The structure of conventional amorphous semiconduc- 
tors like amorphous Si or Ge is well represented by the 
continuous random network model (CRN) introduced by 
ZachariasenlQ] 70 years ago. The CRN model has the sim- 
plicity that each of the atoms should satisfy its local bond- 
ing requirements, and should have minimum strain, charac- 
terized by having a narrow bond angle and bond-length dis- 
tribution. Recently there has been an amorphous graphene 
CRN model proposed 0. Here we develop Kapko et al.'s 
work [ 2], and show that pentagons induce curvature in the 
free standing sheets and analyze the electronic properties. 

2 Models To calculate the band structure, the primi- 
tive cell of graphene was used. To calculate the density of 
states of crystalline graphene, a 800-atom model (800 c- 
g) was made. For amorphous graphene, we used one 800- 
atom (800 a-g) model and two 836-atom models (836 a-gl 
and 836 a-g2). These amorphous graphene models were 
all prepared by a modified Wooten-Weaire- Winer (WWW) 
method0. 



3 Crystalline graphene 

3.1 Band structure When calculating the band struc- 
ture of graphene, tight-binding and ab initio methods are 
two widely used tools. Reich et al. have compared the re- 
sult of tight binding with ab initio[3\. Nowadays, two of 
the widely used ab initio programs are SIESTAQ, using 
pseudopotentials and the Perdew-Zunger parameterization 
of the local-density approximation (LDA), and VASP0, 
with pseudopotentials, plane- wave basis and LDA. 

We computed the eight lowest-energy bands of graphene, 
by using a single-^ (SZ) basis set and a double-^, polarized 
(DZP) basis set with SIESTA, and also VASP. 20 kpoints 
along each special symmetry line were taken for both SZ 
and DZP calculations by SIESTA, and 50 kpoints along 
each line for VASP. The results of SZ and DZP are essen- 
tially identical for occupied bands, and exhibit differences 
for the unoccupied states. Figure 0] shows SIESTA results 
using DZP basis and VASP. We compared our results with 
the other first principle calculations (4] 0. The VASP 
and SIESTA results are in good agreement with published 
results for each code ||4j Q . However, as shown in Figure 
Q] the calculation based on DZP is in good agreement with 
plane-wave pseudopotential calculation (VASP) for the 
four occupied low-energy bands, unlike the unoccupied 
higher energy range, where the agreement is rather poor. 
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Figure 1 Band structure of graphene. a) In SIESTA, the 
solid line represents the result of DZP basis, b) The result 
of VASP is given by dotted-dashed line. 

This discrepancy in the unoccupied bands is interesting, 
and we will explore this further elsewhere. 

Computationally speaking, VASP is more time-consuming 
than SIESTA, particularly for large amorphous graphene 
models we discuss later. For computing total energies and 
forces, e.g. utilizing quantities soley from the occupied 
electronic subspace, SIESTA in SZ approximation is a 
reasonable choice. 

4 Amorphous graphene The three amorphous 
graphene models are prepared by introducing Stone- Wales 
defects into perfect honey comb lattice. The configura- 
tion of 800 a-g model is presented in 0, each atom is 
three-fold. 

4.1 Density of states The electronic density of states 
for the initial flat 800 a-g model is compared to a r point 
density of states for the crystalline 800 c-g model in Figure 
[2 using SIESTA with a SZ basis and the Harris functional. 
From this figure, we observe that the electronic structure 
of 800 a-g model is vastly different from the crystalline 
graphene near the Fermi level, as first reported earlier by 
Kapko et al.[2]. We have constructed additional models 
with periodic boundary conditions and 836 atoms each 
(836 a-gl and 836 a-g2 models). The ring statistics of these 
three models are given in Table Q] and these show some 
small differences. 

Table 1 Ring statistics of 800 a-g, 836 a-gl and 836 a-g2 
models, shown as %. 

Ring Size 800 a-g 836 a-gl 836 a-g2 



5 


33.5 


25 


24 


6 


38 


53 


52 


7 


24 


19 


25 


8 


4.5 


3 






Figure 2 Density of states of 800-atom amorphous and 
crystalline graphene, the Fermi energy is OeV. a) Dot-dash 
line shows the result of 800 a-g model using Harris func- 
tional approximation by SIESTA, c) The density of states 
of 800 c-g model is shown by the dashed line, at the r 
point of Brillouin zone (BZ). The ringing in the crystal is 
due to incomplete BZ sampling. 



4.2 Loss of planar symmetry In all three amor- 
phous graphene models, we introduced small random fluc- 
tuations in the coordinates, in the direction normal to the 
graphene plane, and then relaxed with the Harris func- 
tional and a SZ basis set. Starting with a flat sheet, the 
planar symmetry breaks with curvature above or below 
initial the plane. The final distortion depends on the ini- 
tial conditions. However, a consistent theme emerges of 
pentagons inducing curvature as we describe below. 

As shown in Table [2 Table [3] and Table SI first we 
randomly moved the atoms along normal direction in the 
range of [—Sr, +5r], as shown in the first column of these 
tables; and the results of relaxing by SIESTA in SZ ba- 
sis are shown in the second, third and fourth columns. The 
prime symbol refers to the relaxed model. Taking the 800 
a-g model as an example, the influence of puckering the 
system on the density of states around Fermi level is shown 
in Figure [3j an intuitive view of the fluctuation after relax- 
ing is shown in Figure |4] and [5] when Sr = 0.05A Af- 
ter breaking the planar symmetry by a tiny amount, say 
Sr = 0.05 A, all three models pucker and form the rippled 
or undulated structure as shown in Figure |4] and [5] 

Table 2 the influence of Sr on 800 a-g system relative to 
initial flat model 



Sr (A) 


SF(A) 


Etot/N atom (eV) 


0.01 


0.520 


-0.107 


0.05 


0.525 


-0.107 


0.07 


0.526 


-0.107 
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Table 3 the influence of Sr on 836 a-gl system relative to 
initial flat model 



Sr (A) Sr' (A) E tot /N atom (eV) 



0.01 
0.05 
0.07 



2.53E-3 

1.402 

1.401 



0.0 

-0.102 
-0.102 



Table 4 the influence of Sr on 836 a-g2 system relative to 
initial flat model 



Sr (A) Sr' (A) E tot /N atorn (eV) 



0.01 
0.05 
0.07 



2.72E-3 

1.183 

1.180 



0.0 

-0.090 
-0.090 



0.05 



0.04 



0.03 
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Figure 3 Density of states of the original and relaxed crin- 
kled system, a) The solid line is the density of states of 
original 800-atom amorphous graphene model, b) The den- 
sity of states of crinkled systems are shown as marked in 
the plot, c) The Fermi level is corrected to OeV in the plot, 
as shown in dot- slash line. 




Figure 4 The flat view of the relaxed 836 a-gl system (in 
gray). The blue background illustrates the original 836 a- 
gl model. 




5.809 A 



Figure 5 The side view of the relaxed 800 a-g system. The 
biggest separation along normal direction is marked in the 
plot. 



The radial distribution function g(r) is shown in Fig- 
ure [6] From this plot, the mean bond length of the relaxed 
systems with different initial Sr remain near 1.42 A, the 
change in ring statistics after relaxing is also not signifi- 
cant. And according to Fig HI only one bond broke after re- 
laxation. This implies that the main difference between the 
original amorphous graphene model and the relaxed ones 
is due to these undulations in Figure |U 

To compare with the puckering of the 800 a-g model, 
we also introduced the same planar symmetry breaking 
into a 800 c-g model and relaxed it. As expected, the atoms 
in this crystalline system maintained planar symmetry. 

Also in order to find the relation between the ripples in 
the relaxed systems and the initial random distortion, we 
tested different seeds in random number generator (RNG) 
and also different RNG. The results reveals that the chang- 



ing seeds or or employing different RNG are quantitatively 
small: The maximum mean distortion from the original flat 
plane (Sr') is about 0.545 A and the maximum change in to- 
tal energy is around O.OleV per atom. Figure shows the 
side view of the final configurations by using new and orig- 
inal RNG when Sr = 0.0 5 A, we can tell that the rippled 
regions are similar, except certain regions have formed 
"bucky domes" on opposite sides of the initial plane. 

To further test the relation between the ripples and the 
initial distortion, instead of randomly moving all the atoms 
of 800 a-g, we only distorted the atoms within pentagons 
and compared with atoms not included in pentagons. Fig- 
ure [8] shows the side view of the final configuration of 
relaxed 800 a-g system when the Sr = 0.01 A and only 
the atoms within pentagons were randomly moved . These 
results are similar as the previous test: a) The maximum 
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Figure 6 Radial distribution function of flat and crinkled 
800 a-g system. 



Figure 8 The side view of the final configuration with 
Sr = 0.05 A of only moved atoms within pentagons and 
the original relaxation (distort all atoms), a) The gray balls 
and sticks show the result of moving atoms within pen- 
tagons, b) The blue frames represent the result of original 
distortion. 




Figure 7 The side view of the final configuration by using 
new and original RNG. a) The gray balls and sticks show 
the result of new RNG. b) The blue frames represent the 
result of original RNG. 

change in Sr f is around 0.625A and the maximum change 
in total energy is around 0.02ey per atom, b) No matter 
which atoms were distorted initially, the final puckered re- 
gions involve the same atoms, but possibly puckered in 
the opposite direction relative to the symmetry plane. Fi- 
nally, we note that a 128-atom amorphous graphene model 
made with "melt quenching" [ITU exhibits regions puck- 
ered around pentagons in a similar fashion to what we re- 
port here. 

Different initial symmetry breaking leads to different 
nearly degenerate states after relaxing. However, as stated 
above, no matter how different initial condition is ( or how 
different these degenerated state is), the puckered regions 
are almost the same. It is evident that different rings induce 



these ripples. With this motivation, we searched for regions 
where the height differences of two neighbor atoms are the 
largest and smallest in the model (crinkled and smooth re- 
gions), as shown in Fi gure l9l ITOl and ITT1 In these plots, the 
gray atoms are the configuration of crinkled system, and 
the blue straight lines represent the original model. 

As illustrated in Figure[9]and[T0] the puckered areas are 
pentagon-dense areas. The bonds with most distortion do 
not belong to these pentagons, instead they are within the 
hexagons or heptagons connecting two pentagons. And the 
smooth areas have fewer pentagons than the crinkled ar- 
eas, and most parts of the smooth areas contain hexagons 
and pentagons. These plots imply that the hexagons and 
heptagons alone will not lead to planar symmetry break- 
ing. These ripples formed by pentagons strongly remind 
us of the fullerenes, especially the buckyball (Cqo) which 
only contains pentagons and hexagons. The distance from 
the top to the bottom of the ripples for 800 a-g is around 
5. 809 A as shown in Figure [3 which is comparable to the 
diameter of buckyball, 6.636A. As shown in Figure l9l ITOl 
and Qj] the crinkled regions are all associated with pen- 
tagons. 

5 Conclusion We have shown that random network 
models of amorphous graphene pucker. While the planar 
conformation is locally stable, a lower energy solution is 
obtained that is puckered with local maxima and minima 
in the vicinity of pentagons. The relaxation is performed 
using a density functional calculation of the electronic en- 
ergy. The scale of the puckering is consistent with the cur- 
vature found in buckyball caps with a pentagon surrounded 
by larger rings. 
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(a) (b) 




(c) (d) 



Figure 9 The enlarged plot of crinkled and smooth region 
of 800 a-g model, a) The top view of the crinkled region, 
b) The side view of the crinkled region, c) The top view of 
the smooth region, d) The side view of the smooth region. 




(a) (b) 




(c) (d) 

Figure 1 1 The enlarged plot of crinkled and smooth region 
of 836 a-g2 model, a) The top view of the crinkled region, 
b) The side view of the crinkled region, c) The top view of 
the smooth region, d) The side view of the smooth region. 




(a) (b) 




Figure 1 The enlarged plot of crinkled and smooth region 
of 836 a-gl model, a) The top view of the crinkled region, 
b) The side view of the crinkled region, c) The top view of 
the smooth region, d) The side view of the smooth region. 
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